//==============================================================================
// ReadingFrame.cc
// ----------------------------------------------------------------------------
//
//------------------------------------------------------------------------------
// $Id: ReadingFrame.cc,v 1.1 2000/01/31 21:15:30 kdurbin Exp $
//------------------------------------------------------------------------------
// $Log: ReadingFrame.cc,v $
// Revision 1.1  2000/01/31 21:15:30  kdurbin
// Initial revision
//
//   
//==============================================================================
//.......1.........2.........3.........4.........5.........6.........7.........8
//

#include "ReadingFrame.hh"							

														
// The genetic code
char mCode[4][4][4]={	{{'F','F','L','L'},{'S','S','S','S'},{'Y','Y','.','.'},{'C','C','.','W'}},			
					 	{{'L','L','L','L'},{'P','P','P','P'},{'H','H','Q','Q'},{'R','R','R','R'}},
					  	{{'I','I','I','M'},{'T','T','T','T'},{'N','N','K','K'},{'S','S','R','R'}},
					  	{{'V','V','V','V'},{'A','A','A','A'},{'D','D','E','E'},{'G','G','G','G'}}};


// 022, for example, corresponds to {U/T}AA which is a stop

//---------------------------------------------------------------------------
// * Idx
//---------------------------------------------------------------------------
// Returns a numerical index corresponding to the specified letter. 
// U or T  = 0
// C       = 1
// A       = 2
// G       = 3
// unknown = 4
//
int
ReadingFrame::Idx(char c)
{
	switch(c){
	case 'U':
	case 'T':
	case 'u':
	case 't':
		return(0);
		break;
	case 'c':
	case 'C':
		return(1);
		break;
	case 'a':
	case 'A':
		return(2);
		break;
	case 'g':
	case 'G':
		return(3);
		break;
	default:
		return(4);
		break;
	}
}

//---------------------------------------------------------------------------
// * Decode
//---------------------------------------------------------------------------
//
// 
char
ReadingFrame::Decode(char a0,char a1,char a2)
{
	int i0=Idx(a0),i1=Idx(a1),i2=Idx(a2);
	if ((i0==4) || (i1==4) || (i2 == 4)){
		return('?');
	}
	return(mCode[i0][i1][i2]);
}

//---------------------------------------------------------------------------
// * DecodeSeq
//---------------------------------------------------------------------------
// 
void
ReadingFrame::DecodeSeq(int start,int frame,
						Sequence &DNAIn,
						Sequence &ProteinOut)
{
	int DIdx;
	for(DIdx = start+frame;DIdx < DNAIn.size();DIdx+=3){
		char a0 = DNAIn[DIdx];
		char a1 = DNAIn[DIdx+1];
		char a2 = DNAIn[DIdx+2];
		ProteinOut.push_back(Decode(a0,a1,a2));
	}
}

//---------------------------------------------------------------------------
// * IsStop
//---------------------------------------------------------------------------
// true if codon is one of:
// UAA UAG UGA
// TAA TAG TGA
//
// (U/T)A(A/G)  (U/T)GA
//
bool
ReadingFrame::IsStop(char c1,char c2,char c3)
{	
	if(Decode(c1,c2,c3) == '.') return(true);
	else return(false);
}

//---------------------------------------------------------------------------
// * IsStart
//---------------------------------------------------------------------------
// true if codon is MET (ATG)
//
bool
ReadingFrame::IsStart(char c1,char c2,char c3)
{
	if(Decode(c1,c2,c3) == 'M') return(true);
	else return(false);
}


//---------------------------------------------------------------------------
// * FindFwdOpenReadingFrames
//---------------------------------------------------------------------------
// Find the open reading frames in the forward direction.  Looks for MET-STOP
// intervals.  
//
void
ReadingFrame::FindOpenReadingFrames(Sequence &DNAIn,
									ORFVec &ORFOut,
									int start,int end,
									int minSize,bool IsFwd)
{
	int lastStop = 0;
	OpenReadingFrame temp;
	
	// Forward frames are numbered 0,1,2, reverse 3,4,5.  
	int sense = 0;
	if (!IsFwd){
		sense = 3;
	}

	// Try all three reading frames
	for(int frame = 0;frame < 3;frame++){	
		// Indicate that we are waiting to see a MET
		temp.start = -1;
		// Go through the string, recording ORFs	
		for(int dIdx = start+frame;dIdx < end;dIdx+=3){
			// Record FIRST MET we find only.  Then wait for 
			// a STOP.  
			if (IsStart(DNAIn[dIdx],DNAIn[dIdx+1],DNAIn[dIdx+2])&&
				(temp.start == -1)){
				temp.start = dIdx;
			}else{	
				if (IsStop(DNAIn[dIdx],DNAIn[dIdx+1],DNAIn[dIdx+2])){
					// If we've seen another MET
					if (temp.start != -1){
						// Length is after MET to before . 
						temp.length = dIdx - temp.start-3; 
						temp.frame = frame+sense;
						if (temp.length >= minSize){
							ORFOut.push_back(temp);
						}
						temp.start = -1;
					}
				}
			}
		}
	}
}


//---------------------------------------------------------------------------
// * FindOpenReadingFrames2
//---------------------------------------------------------------------------
// Find the open reading frames in the forward direction.  This version
// looks for stop-stop intervals. 
//
void
ReadingFrame::FindOpenReadingFrames2(Sequence &DNAIn,
									 ORFVec &ORFOut,
									 int minSize,bool IsFwd)
{
	int lastStop = 0;
	OpenReadingFrame temp;
	
	// Forward frames are numbered 0,1,2, reverse 3,4,5.  
	int sense = 0;
	if (!IsFwd){
		sense = 3;
	}

	// Try all three reading frames
	for(int frame = 0;frame < 3;frame++){
		// Go through the string, recording ORFs
		for(int dIdx = frame;dIdx < DNAIn.size();dIdx+=3){
			if (IsStop(DNAIn[dIdx],DNAIn[dIdx+1],DNAIn[dIdx+2])){
				// Special case when lastStop is beginning of seq. 
				if (lastStop > frame){
					temp.start = lastStop+4;
					temp.length = dIdx - temp.start+1; // After . to Before . 
				}else{
					temp.start = lastStop;
					temp.length = dIdx - temp.start; // beginning to Before .
				}
				temp.frame = frame+sense;
				if (temp.length >= minSize){
					ORFOut.push_back(temp);
				}
				lastStop = dIdx;
			}
		}
	}
}


//---------------------------------------------------------------------------
// * FindFwdOpenReadingFrames
//---------------------------------------------------------------------------
// Find the open reading frames in the forward direction.  
//
void
ReadingFrame::FindFwdOpenReadingFrames(Sequence &DNAIn,
										ORFVec &ORFOut,
										int minSize)
{
	FindOpenReadingFrames(DNAIn,ORFOut,0,DNAIn.size(),minSize,true);
}


//---------------------------------------------------------------------------
// * FindAllOpenReadingFrames
//---------------------------------------------------------------------------
// Find open reading frames in all six senses. 
//
void
ReadingFrame::FindAllOpenReadingFrames(Sequence &DNAIn,
										ORFVec &ORFOut,
										int minSize)
{
	// Find ORFs in 3 forward frames
	FindOpenReadingFrames(DNAIn,ORFOut,0,DNAIn.size(),minSize,true);
	
	// Find ORFs in 3 reverse frames
	DNAIn.ReverseComplement();
	FindOpenReadingFrames(DNAIn,ORFOut,0,DNAIn.size(),minSize,false);
	
	// Restore sequence to original state. 
	DNAIn.ReverseComplement();
}



void
ReadingFrame::FindAllStarts(Sequence &DNAIn,ORFVec &ORFOut){
  FindStarts(DNAIn,ORFOut,true);
  
  // Find reverse complement....
  DNAIn.ReverseComplement();
  FindStarts(DNAIn,ORFOut,false);
  DNAIn.ReverseComplement();
}

void
ReadingFrame::FindAllStops(Sequence &DNAIn,ORFVec &ORFOut){
  FindStops(DNAIn,ORFOut,true);
  
  // Find reverse complement....
  DNAIn.ReverseComplement();
  FindStops(DNAIn,ORFOut,false);
  DNAIn.ReverseComplement();
}


//---------------------------------------------------------------------------
// * FindStarts
//---------------------------------------------------------------------------
// Finds all MET codons, in any of the three fwd reading frames. 
//
void
ReadingFrame::FindStarts(Sequence &DNAIn,ORFVec &ORFOut,bool IsFwd)
{
  OpenReadingFrame temp;
	
  // Forward frames are numbered 0,1,2, reverse 3,4,5.  
  int sense = 0;
  if (!IsFwd){
	 sense = 3;
  }

  // Try all three reading frames
  for(int frame = 0;frame < 3;frame++){	
	 // Go through the string, recording ORFs	
	 for(int dIdx = frame;dIdx < DNAIn.size();dIdx+=3){
		// Record each MET we find. 
		if (IsStart(DNAIn[dIdx],DNAIn[dIdx+1],DNAIn[dIdx+2])){
				temp.start = dIdx;
				temp.length = 3;
				temp.frame = frame+sense;
				ORFOut.push_back(temp);
		}
	 }
  }
}


//---------------------------------------------------------------------------
// * FindStops
//---------------------------------------------------------------------------
// Finds all STOP codons, in any of the three fwd reading frames. 
//
void
ReadingFrame::FindStops(Sequence &DNAIn,ORFVec &ORFOut,bool IsFwd)
{
  OpenReadingFrame temp;
	
  // Forward frames are numbered 0,1,2, reverse 3,4,5.  
  int sense = 0;
  if (!IsFwd){
	 sense = 3;
  }

  // Try all three reading frames
  for(int frame = 0;frame < 3;frame++){	
	 // Go through the string, recording ORFs	
	 for(int dIdx = frame;dIdx < DNAIn.size();dIdx+=3){
		// Record each MET we find. 
		if (IsStop(DNAIn[dIdx],DNAIn[dIdx+1],DNAIn[dIdx+2])){
				temp.start = dIdx;
				temp.length = 3;
				temp.frame = frame+sense;
				ORFOut.push_back(temp);
		}
	 }
  }
}

